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Abstract Using Path Integral Monte Carlo and the Max- 
imum Entropy method, we calculate the dynamic struc- 
ture factor of solid ^He in the bcc phase at a finite 
temperature of T = 1.6 K and a molar volume of 21.5 
cm^. From the single phonon dynamic structure factor, 
we obtain both the longitudinal and transverse phonon 
branches along the main crystalline directions, [001], [Oil] 
and [111]. Our results are compared with other theoret- 
ical predictions and available experimental data. 

Keywords quantum solid • Path Integral Monte Carlo • 
Maximum Entropy method ■ phonon spectrum • solid 
helium 
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1 Introduction 

The two isotopes of solid helium, ^He and '^He, are the 
most prominent examples of quantum solids |21E1 ■ These 
isotopes display highly anharmonic dynamics: their atoms 
are strongly correlated, loosely bound and make large 
excursions to the nearest neighbor sites. The effects are 
especially significant in the bcc phase of solid helium. 
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The dynamics of solid '^He, in particular the phonon 
spectrum, has been studied over years in inelastic neu- 
tron scattering experiments 01151017] • In contrast to ^He, 
there are no experimental measurements of phonon spec- 
tra of '^He (in the bcc phase) due to the large neutron 
absorption cross-section of this isotope. OIHl Recently, 
however, the use of inelastic X-ray scattering |^10| offers 
the possibility of experimentally determining the phonon 
branches of ^He. In addition, new neutron scattering ex- 
periments are being planned, |11II12| which will be able 
to acquire data for a sufficiently long time to observe in- 
elastic scattering. The present study has been motivated 
by these developments. 

Strong anharmonic effects and atomic correlations 
make the theoretical calculation of phonon spectra of 
•^He very difficult. Several theoretical methods have been 
developed to treat this problem: Clyde and Khanna de- 
vised a self-consistent phonon t-matrix formalism |13|. 
Horner applied a many-body perturbation technique [141 
115) , and Koehler and Werthamer employed a variational 
approach. These theoretical calculations are based on a 
variational perturbative theory and implemented at zero- 
temperature. (21 The predictions of these models are in 
quantitative disagreement with each other, due to differ- 
ences in the effective potentials used to evaluate phonon 
interactions and details of the numerical methods em- 
ployed, m 

As a complementary approach to these methods we 
decided to study the excitations in bcc solid helium ^He, 
by performing Quantum Monte Carlo simulations at a 
finite temperature. We use Path Integral Monte Carlo 
(PIMC) XLn which is a non-perturbative numerical method, 
that allows, in principle, simulations of quantum systems 
without any assumptions beyond the Schrodinger equa- 
tion. The two body interatomic He-He potential ^S] is 
the only input for the PIMC simulations. In our study 
the Universal Path Integral code of Ceperley^^ '^^s 
adapted to calculate the phonon branches at finite tem- 
perature. The PIMC method, in conjunction with Max- 
imum Entropy (MaxEnt) techniques, (19^ was employed 
to calculate the phonon spectra of liquid [211 and solid [21] 
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helium ^He) at finite temperatures. The calculated phonon 
spectra are in good agreement with the spectra measured 
by inelastic neutron scattering experiments. 

The results of our study include the numerical calcu- 
lation of all the phonon branches of bcc ^He at a molar 
volume 21.5 cm'^. Details of our simulations are described 
in Sec. 2, the results of the calculations are presented in 
Sec. 3. 



2 Method 

The PIMC method used in our simulations is based on 
the formulation of quantum mechanics in terms of path 
integrals. It has been described in detail by Ceperley 
|17) . The method involves mapping of the quantum sys- 
tem of particles onto a classical model of interacting 
"ring polymers" , whose elements, "beads" or "time-slices" , 
are connected by "springs". The method provides a di- 
rect statistical evaluation of quantum canonical aver- 
ages. In addition to static properties of the system, dy- 
namical properties can be also extracted from PIMC 
simulations. '17' 

The object of this study is the phonon spectrum, 
which can be extracted from the dynamic structure fac- 
tor, ^(q, uj). The definition of S{q, uj) in terms of density 
fluctuations is 
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where ?iq and huj are the momentum and energy (we 
take ?i = 1), pq is the Fourier transform of the density of 
the solid, and n is the number density, ^(q, to) is usually 
expressed in terms of phonons, by writing S{q,Lu) as a 
sum of terms involving the excitation of a single phonon, 
S'i(q, w), a pair of phonons, S'2(q, ti-") and higher order 
terms which also include interference between different 
terms. pill4| In most of our simulations we calculated 
the Si (q, oj) term. Some calculations of S'(q, lu) were also 
performed, and will be discussed below. Additional de- 
tails of the method used to calculate both 5*1 (q,a;) and 
S{q,uj) and extract phonon spectra are described in our 
paper. |22 

In the simulations we used samples containing be- 
tween 128 and 432 atoms. Each atom was represented 
by a "ring" polymer with 64 time slices. S'(q, to) was cal- 
culated for values of q between 0.14 and 1 in relative 
lattice units (r.l.u.= 2Tr/a, where a is the lattice param- 
eter). The number density was set to p = 0.02801 (1/A^) 
and the temperature to T=1.6 K. This particular tem- 
perature was chosen to compare the phonon spectrum of 
■^He with those of '*He calculated at the same tempera- 
ture A perfect bcc lattice was prepared for the ini- 
tial configuration. The effects of Fermi statistics are not 
taken into account in our simulations, which is a reason- 
able approximation for the solid phase at T=1.6 K. 



Statistical errors were estimated by running the PIMC 
simulations at least 10 times, with different initial condi- 
tions in each case. After each run, S{q,u>) (5'i(q, w)) was 
extracted using the MaxEnt method HH]. The pho non 
energy for a given q was then calculated by averaging the 
positions of the peak of S'(q, uj) over the set of the sim- 
ulation runs. The error bars shown in the figures below 
represent the standard deviation. We collected at least 
10000 data points in each simulation run and re-blocked 
the data |2Ij in blocks of 100-200 points |2Ij. Each sim- 
ulation run took about two weeks of 12 Pentium III PCs 
running in parallel. 

A typical example of the calculated dynamic struc- 
ture factor 5'i(q, along the [Oil] direction is shown 
in Fig^ The figure shows both the longitudinal and 
transverse phonons. As seen in Fignthe two transverse 
phonons have narrow and almost symmetric line shapes, 
while the longitudinal phonon is quite broad and asym- 
metric. This asymmetry exhibits itself in the form of a 
relatively steep rise on the low-frequency side and a more 
gradual decrease on the high-frequency side. Similar form 
of the phonon spectral function was obtained by Koehler 
and Wcrthamer .16.. 




Fig. 1 Longitudinal component (solid line) and transverse 
components (Ti - dotted line, T2 - dashed line) of the dy- 
namic structure factor, Si{q,ij), for q — 0.4 r.l.u. along the 
[Oil] direction. 



3 Results 

The calculated dispersion relations of longitudinal and 
transverse phonons, along the main crystal directions 
([001], [111] and [Oil]), are shown in Figs. (The 
numerical values are given in the Appendix). 

For comparison, we plot the phonon branches calcu- 
lated at T=0 K by Clyde and Khanna by Koehler 
and Wcrthamer 16^ (at molar volume V=21.5 cni^^) and 
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0.4 0.6 
q (r.l.u) 

Fig. 2 Calculated dispersion relation of the L[001] phonon 
branch (triangles) using 5*1 (q, a;). The error bars represent 
statistical uncertainty. For refe renc e, phonon frequencies cal- 
culated by Clyde and Khanna 1131 (solid lines), Koehler and 
Werthamer |1d| (dotted lines) and Horner 1 5_ (dashed line) 
are shown. The phonon spectrum at small q (circles) is es- 
timated using elastic constants measured by Grey wall ,22 . 
The scaled single-phonon branch of ''He (dot-dashed line) is 
taken from |21|. 




0.4 0.6 
q (r.Lu.) 

Fig. 3 Calculated dispersion relation of the T[001] phonon 
branch (triangles) using Si{q,uj). The error bars represent 
statistical uncertainty. For reference, phonon frequencies cal- 
culated by Clyde and Khanna 1131 (solid lines), Koehler and 
Werthamer (1d| (dotted lines) and Horner I_5 (dashed line) 
are shown. The phonon spectrum at small q (circles) is es- 
timated using elastic constants measured by Grey wall 1221. 
The scaled single-phonon branch of ''He (dot-dashed line) is 
taken from 1211. 




0.2 0.3 
q (r.l.u.) 

Fig. 4 Calculated dispersion relation of the L[011] phonon 
branch (triangles) using S'i(q, tj). The error bars represent 
statistical uncertainty. For refe renc e, phonon frequencies cal- 
culated by Clyd e and Khanna tl3| (solid lines), Koehler and 
Werthamer |lti| (dotted lines) and Horner [T5| (dashed line) 
are shown. The scaled single-phonon branch of *He 
dashed line) is taken from 21 . 
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0.2 0.3 
q (r.l.u.) 

Fig. 5 Calculated dispersion relations of transverse phonon 
branches along [Oil] using 5'i(q, aj). Calculated values are 
shown for the Ti branch (triangles down) and T2 branch (tri- 
angles up). The error bars represent statistical uncertainty. 
For reference, phonon frequencies calculated by Clyde and 
Khanna 13 (solid lines), Koehler and Werthamer (dotted 
lines) and Horner 15: (dashed line) are shown. The phonon 
spectra at small q (circles) ar e es timated using elastic con- 
stants measured by Creywall |22| . The scaled single-phonon 
branches of ''He (dot-dashed lines) is taken from |21|. 



Horner ^S] (at molar volume V=24 cm^). Koehler and 
Werthamer suggested that due to the asymmetry of the 
dynamic structure factor, the phonon frequency can be 
calculated either at the position of the maximum of ^(q, co) 
or as the mean of the two half-maxima jl6) . Since we 
calculated the phonon spectra using the position of the 



maximum of S'(q, lu) only these phonon frequencies were 
taken from Koehler and Werthamer |TB| . 

In addition, we evaluated phonon frequencies in the 
long wavelength limit (q — + 0) using the elastic moduli 
of bcc ^He measured by Greywalll55l at molar volume V 
= 21.6 cm'^. We also made a fit to ^He phonon branches 
calculated by the PIMC method scaling them by 
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0.2 0.4 0.6 0.8 1 
q (r.l.u.) 

Fig. 6 Calculated dispersion relation of the L[lll] phonon 
branch (triangles) using 5'i(q, oj). The error bars represent 
statistical uncertainty. The error bars represent statistical 
uncertainty. For refe renc e, phonon frequencies calculated by 
Clyde and Khanna [131 (solid lines) along with frequencies 
calculated by Koehler and Werthamer (16. (dotted lines) are 
shown. The scaled sing le-phonon branch of ^He (dot-dashed 
line) is taken from '21 . 




0.2 0.4 0.6 0.8 1 
q (r.l.u.) 

Fig. 7 Calculated dispersion relation of the T[lll] phonon 
branch (triangles) using 5*1 (q, a;). The error bars represent 
statistical uncertainty. The error bars represent statistical 
uncertainty. For refe renc e, phonon frequencies calculated by 
Clyde and Khanna 1131 (solid lines) al ong with frequencies 
calculated by Koehler and Werthamer |16|(' dotted lines) are 
shown. The scaled sing le-phonon branch of *He (dot-dashed 
line) is taken from |21|. 

square root of m('*He)/ni(^He), the mass ratio of the 
isotopes, and plot in Figs.|3-[3 

Surprisingly, we found that the results of ours sim- 
ulations show good agreement with the theoretical cal- 
culation of Glyde and Khanne (except the phonon 
branch Ti [Oil] in Fig. In addition, as can be seen 
from Figs. |2]-|7| the calculated phonon branches of "^He 
coincide very closely with the scaled branches of "^He 
(being practically on top of each other for longitudinal 



phonons). We conclude that, according to our simula- 
tions, the phonon frequencies of these two isotopes scale 
as an inverse ratio of the square root of their atomic 
masses. 



4 Conclusions 

In conclusion, we calculated the dynamic structure fac- 
tor for solid ^He in the bcc phase of using PIMC simu- 
lations and the MaxEnt method. PIMC was used to cal- 
culate the intermediate scattering function in the imag- 
inary time from which the dynamic structure factor was 
inferred with the MaxEnt method. We extracted both 
the longitudinal and transverse phonon branches from 
the one-phonon dynamic structure factor. Good agree- 
ment between the calculated phonon spectra and the 
theoretical prediction of Glyde and Khanne ^31 has been 
obtained. The results of our simulations show that the 
phonon frequencies of these two isotopes scale as an in- 
verse ratio of the square root of their atomic masses. 
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5 Appendix: Phonon energies 

The energies of phonons calculated from S'i(q, cj) in the 
bee phase of "^He (molar volume 21 cm'^) and '^He (molar 
volume 21.5 cm3) at T = 1.6 K are listed in Tabs.^H 
The phonon energy, uj, is in units of Kelvin, and the 
reciprocal lattice vector, q, is in relative lattice units 
(r.l.u.= 27r/a, where a = 4.1486A is the lattice parame- 
ter). 



Table 1 Calculated phonon energies, of the L[001], 

T[001], L[lll] and T[lll] phonon branches of bcc *He ob- 
tained by using 5'i(q,tj). 



Table 3 Calculated phonon energies a;(K) of the L[001], 
T[001], L[lll] and T[lll] phonon branches of bcc ^He ob- 
tained by using 5i(q,a;). 



q (r.I.u) 


L[001] 


T[001] 


L[lll] 


T[lll] 


0.20 


12.6 


± 


0.5 


8.7 ± 0.5 


19.2 


± 


2.9 


9.7 ± 0.4 


0.25 


15.8 


± 


0.5 


10.5 ± 0.5 


24.5 


± 


3.1 


11.7 ± 0.6 


0.33 


17.0 


± 


1.6 


11.9 ± 1.1 


23.9 


± 


3.9 


14.3 ± 1.7 


0.40 


22.0 


± 


1.2 


15.9 ± 0.5 


25.8 


± 


4.7 


18.4 ± 0.8 


0.50 


23.6 


± 


0.6 


18.9 ± 0.9 


21.3 


± 


1.1 


21.9 ± 1.5 


0.60 


22.1 


± 


1.5 


21.3 ± 1.0 


13.5 


± 


1.5 


24.7 ± 2.1 


0.67 


23.5 


± 


2.4 


22.0 ± 1.7 


11.9 


± 


0.7 


25.8 ± 3.5 


0.75 


23.5 


± 


1.6 


24.5 ± 1.1 


15.8 


± 


3.1 


24.3 ± 2.1 


0.80 


24.9 


± 


0.9 


25.5 ± 2.0 


15.0 


± 


3.6 


24.7 ± 1.9 


0.83 


24.6 


± 


1.4 


25.2 ± 1.2 


19.4 


± 


4.1 


25.2 ± 2.1 


1.00 


26.6 


± 


0.8 


25.4 ± 1.3 


25.7 


± 


0.7 


25.7 ± 2.7 



Table 4 Calculated phonon energies tt'(K) of the L[011], 
T2[011] and Ti[011] phonon branches of ^He obtained by us- 
ing Si{ci,uj). 



q (r.I.u) 


L[001] 


T[001] 


L[lll] 


T[lll] 


0.14 


7.1 ± 1.9 


6.1 ± 1.0 


16.2 ± 1.5 


6.1 ± 1.9 


0.20 


9.1 ± 2.5 


8.1 ± 1.5 


20.2 ± 2.8 


8.1 ± 0.0 


0.25 


12.1 ± 1.5 


9.1 ± 1.8 


21.2 ± 2.2 


10.1 ± 1.8 


0.29 


13.1 ± 1.6 


10.1 ± 1.5 


22.2 ± 4.0 


10.1 ± 1.5 


0.33 


15.5 ± 1.4 


12.1 ± 1.0 


22.2 ± 1.2 


13.1 ± 1.8 


0.40 


16.2 ± 1.6 


13.1 ± 1.1 


23.2 ± 2.7 


15.2 ± 2.1 


0.43 


18.2 ± 1.5 


15.2 ± 1.8 


23.2 ± 3.4 


16.2 ± 3.0 


0.50 


19.2 ± 1.1 


16.2 ± 1.2 


18.2 ± 1.9 


18.2 ± 2.2 


0.57 


18.2 ± 1.4 


18.2 ± 1.5 


14.2 ± 4.1 


20.2 ± 2.2 


0.60 


20.2 ± 1.7 


19.2 ± 0.9 


11.1 ± 2.0 


20.2 ± 2.8 


0.67 


21.2 ± 1.7 


19.2 ± 1.8 


9.1 ± 2.0 


22.2 ± 3.4 


0.71 


21.2 ± 1.9 


20.2 ± 1.9 


10.1 ± 1.7 


22.2 ± 4.2 


0.75 


21.2 ± 2.2 


20.2 ± 1.9 


13.1 ± 1.7 


21.2 ± 3.6 


0.80 


21.2 ± 2.2 


21.2 ± 2.4 


15.2 ± 1.7 


22.2 ± 2.3 


0.83 


22.4 ± 1.7 


22.5 ± 2.3 


18.2 ± 1.6 


24.2 ± 4.1 


0.86 


21.2 ± 2.6 


22.2 ± 2.5 


22.2 ± 1.8 


23.2 ± 3.1 



q (r.I.u) 


L[011] 


Ti[011] 


T2[011] 


0.17 


18.5 ± 1.7 


1.8 ± 1.5 


10.0 ± 0.7 


0.20 


20.9 ± 0.5 


2.8 ± 3.8 


11.2 ± 0.3 


0.25 


24.2 ± 1.0 


6.1 ± 2.8 


13.1 ± 0.6 


0.33 


30.1 ± 1.6 


7.1 ± 0.5 


15.8 ± 0.7 


0.40 


32.0 ± 2.5 


7.1 ± 0.0 


17.2 ± 0.4 


0.50 


32.4 ± 2.0 


8.1 ± 2.5 


17.6 ± 1.1 



Table 2 Calculated phonon energies uj(K) of the L[011], 
T2[011] and Ti[011] phonon branches of bcc ^He obtained by 
using S'i(q,a;). 



q (r.I.u) 


L[011] 


Ti[011] 


T2[011] 


0.20 


16.2 ± 2.0 


5.1 ± 1.2 


10.1 ± 1.3 


0.25 


33.3 ± 2.3 


5.1 ± 1.7 


12.1 ± 1.3 


0.29 


21.2 ± 2.7 


6.0 ± 1.7 


12.1 ± 1.0 


0.33 


26.3 ± 1.5 


6.05 ± 1.1 


14.2 ± 2.0 


0.40 


29.3 ± 1.9 


6.11 ± 1.0 


15.2 ± 1.5 


0.43 


26.3 ± 2.8 


6.1 ± 1.6 


16.2 ± 2.5 


0.50 


31.3 ± 2.5 


6.12 ± 2.4 


15.2 ± 1.8 



